/*@HEADER
// ***********************************************************************
//
//       Ifpack2: Templated Object-Oriented Algebraic Preconditioner Package
//                 Copyright (2009) Sandia Corporation
//
// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
// license for use of this work by or on behalf of the U.S. Government.
//
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are
// met:
//
// 1. Redistributions of source code must retain the above copyright
// notice, this list of conditions and the following disclaimer.
//
// 2. Redistributions in binary form must reproduce the above copyright
// notice, this list of conditions and the following disclaimer in the
// documentation and/or other materials provided with the distribution.
//
// 3. Neither the name of the Corporation nor the names of the
// contributors may be used to endorse or promote products derived from
// this software without specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
//
// Questions? Contact Michael A. Heroux (maherou@sandia.gov)
//
// ***********************************************************************
//@HEADER
*/

#include "Teuchos_CommandLineProcessor.hpp"
#include "Teuchos_StandardCatchMacros.hpp"
#include "Teuchos_VerboseObject.hpp"
#include "Teuchos_ParameterList.hpp"
#include "Teuchos_ParameterXMLFileReader.hpp"
#include "Teuchos_Time.hpp"
#include "Teuchos_StackedTimer.hpp"
#include "Teuchos_Comm.hpp"

#include "Ifpack2_Parameters.hpp"

#include "Tpetra_Core.hpp"
#include "Tpetra_Map.hpp"
#include "Tpetra_BlockCrsMatrix_Helpers_decl.hpp"
//#include "belos_header.pch"
//#include "belos_inc.pch"
#include "build_problem.hpp"
#include "build_solver.hpp"
#include <cmath>
//Teuchos::RCP<Belos::SolverManager<Scalar,TMV,TOP> >
//build_solver(Teuchos::ParameterList& test_params,
//             Teuchos::RCP<Belos::LinearProblem<Scalar,TMV,TOP> > problem);

#define RESET   "\033[0m"
#define BLACK   "\033[30m"      /* Black */
#define RED     "\033[31m"      /* Red */
#define GREEN   "\033[32m"      /* Green */
#define YELLOW  "\033[33m"      /* Yellow */
#define BLUE    "\033[34m"      /* Blue */
#define MAGENTA "\033[35m"      /* Magenta */
#define CYAN    "\033[36m"      /* Cyan */
#define WHITE   "\033[37m"      /* White */


void
process_command_line (bool& printedHelp,
                      std::string& xml_file,
                      bool& useStackedTimer,
                      std::string& problem_name,
                      int argc,
                      char*argv[])
{
  Teuchos::CommandLineProcessor cmdp(false,true);
  cmdp.setOption("xml_file", &xml_file, "XML Parameters file");
  cmdp.setOption("with_stacked_timer", "without_stacked_timer", &useStackedTimer,
      "Whether to run with a StackedTimer and print the timer tree at the end");
  cmdp.setOption("problem_name", &problem_name, "Problem name for Watchr plot");

  const auto result = cmdp.parse (argc, argv);

  // mfh 21 Apr 2016: By ignoring options that this executable doesn't
  // recognize, we can pass them through to (e.g.,) Kokkos.

  if (result == Teuchos::CommandLineProcessor::PARSE_HELP_PRINTED) {
    printedHelp = true; // not an error to ask for help
  }
  else if (result == Teuchos::CommandLineProcessor::PARSE_ERROR) {
    throw std::runtime_error ("Error parsing command-line.");
  }
}
int main (int argc, char* argv[])
{
  typedef double Scalar;
  typedef Teuchos::ScalarTraits<Scalar>        STS;
  typedef STS::magnitudeType                   Magnitude;
  typedef Tpetra::Map<>::node_type             Node;
  typedef Tpetra::Map<>::local_ordinal_type    LO;
  typedef Tpetra::Map<>::global_ordinal_type   GO;
  typedef Tpetra::MultiVector<Scalar,LO,GO>    TMV;
  typedef Tpetra::Operator<Scalar,LO,GO>       TOP;
  typedef Tpetra::CrsMatrix<Scalar,LO,GO,Node> crs_matrix_type;
  typedef Tpetra::BlockCrsMatrix<Scalar,LO,GO,Node> block_crs_matrix_type;
  typedef Belos::LinearProblem<Scalar,TMV,TOP> BLinProb;
  typedef Belos::SolverManager<Scalar,TMV,TOP> BSolverMgr;
  using Teuchos::RCP;
  using Teuchos::rcp;
  using Teuchos::StackedTimer;

  Tpetra::ScopeGuard tpetraScope (&argc, &argv);

  bool success = true;

  RCP<Teuchos::FancyOStream>
    out = Teuchos::VerboseObjectBase::getDefaultOStream();

  try {


    RCP<const Teuchos::Comm<int> > comm =
      Tpetra::getDefaultComm();


    //Just get one parameter from the command-line: the name of an xml file
    //to get parameters from.

    std::string xml_file("poly_test.xml");
    std::string problem_name("Belos-Ifpack2 Setup-Solve");
    bool useStackedTimer = false;
    {
      bool printedHelp = false;
      process_command_line (printedHelp, xml_file, useStackedTimer, problem_name, argc, argv);
      if (printedHelp) {
        return EXIT_SUCCESS;
      }
    }

    //Read the contents of the xml file into a ParameterList. That parameter list
    //should specify a matrix-file and optionally which Belos solver to use, and
    //which Ifpack2 preconditioner to use, etc. If there are sublists of parameters
    //for Belos and Ifpack2, those will be passed to the respective destinations
    //from within the build_problem and build_solver functions.

    //*out << "Every proc reading parameters from xml_file: "
    //     << xml_file << std::endl;
    Teuchos::ParameterList init_params =
      Teuchos::ParameterXMLFileReader(xml_file).getParameters();
    
    RCP<StackedTimer> stackedTimer;

    //The build_problem function is located in build_problem.hpp.
    //Note that build_problem calls build_precond and sets a preconditioner on the
    //linear-problem, if a preconditioner is specified.
    int mm_start,mm_end; 
    int prec_start,prec_end, iter_num; 
    bool test_new(false);
    Ifpack2::getParameter(init_params, "mm_start", mm_start);
    Ifpack2::getParameter(init_params, "mm_end", mm_end);
    Ifpack2::getParameter(init_params, "prec_start", prec_start);
    Ifpack2::getParameter(init_params, "prec_end", prec_end);
    Ifpack2::getParameter(init_params, "test_new", test_new);
    Ifpack2::getParameter(init_params, "iter_num", iter_num);
    
    std::set<int> skipIndices = {7, 9, 11, 12, 16, 18, 20, 22, 25, 26, 28, 29, 32, 33, 35,40,41,42,43};
    std::set<int> skipIndices1 ={1,6,7,21,22};//{7, 9, 11, 12, 16, 18, 20, 22, 25, 26, 28, 29, 32, 33, 35,40,41,42,43};
    if(test_new) skipIndices = skipIndices1;
    for (int n = mm_start; n <= mm_end; ++n) {
        if (skipIndices.find(n) != skipIndices.end()) {
            continue; // 跳过该次循环
        }
    //if(n==25) n=n+2;
    for (int iter=1; iter <= iter_num; ++iter) {
    for (int eig=1; eig <= 1; ++eig) {//=1
    for (int prec = prec_start; prec <= prec_end; ++prec) {//<=3
        
      Teuchos::ParameterList test_params =
        Teuchos::ParameterXMLFileReader(xml_file).getParameters();
        Teuchos::Time timer("total");
        timer.start();
        std::stringstream ss;
        if(test_new)     
          ss << "mat_mtx_test/mat_" << n << ".mtx";
        else
          ss << "mat_mtx/mat_" << n << ".mtx";
        std::string filename = ss.str();
        test_params.set("mm_file", filename);

        std::string mm_file;
        Ifpack2::getParameter(test_params, "mm_file", mm_file);
        //if(iter==2 && prec<3)
        //    std::cout <<CYAN<< "============== run again with different order =============="<<RESET<< std::endl;
        //std::cout << "mm_file value for mat_" << n << ".mtx: " << mm_file << std::endl;

        if (prec == 1) {
            std::cout << mm_file <<" ==============================="<< std::endl;
            //std::cout <<YELLOW<< "+++++++++++ Neumann:ritz_compute_power +++++++++++"<<RESET<< std::endl;
            if(iter==2){
            std::cout <<CYAN<< "----------- Neumann:ritz_updateSol_ritz -----------"<<RESET<< std::endl;
            }else{
            std::cout <<YELLOW<< "----------- Neumann:ritz_updateSol_ritz -----------"<<RESET<< std::endl;
            }
            //    test_params.set("solver_type", "GmresPolySolMgr");
            //    test_params.set("Preconditioner Side", "None");
            test_params.set("solver_type", "Block Gmres");
            test_params.set("Ifpack2::Preconditioner", "NEUMANN");
            test_params.set("Preconditioner Side", "Right");

            if (!test_params.sublist("Ifpack2").isParameter("neumann: degree") || iter==2) {
              test_params.sublist("Ifpack2").set("neumann: degree", 2);
            }
            //test_params.sublist("Ifpack2").set("neumann: degree", 1);
            test_params.sublist("Ifpack2").set("neumann: eigenvalue max iterations", 10);
            //test_params.sublist("Ifpack2").set("eigen-analysis: type", "power method");
            if (!test_params.sublist("Ifpack2").isParameter("eigen-analysis: type")) {
              test_params.sublist("Ifpack2").set("eigen-analysis: type", "ritz");
            }
            
            if (!test_params.sublist("Ifpack2").isParameter("neumann: algorithm")) {
              test_params.sublist("Ifpack2").set("neumann: algorithm", "neumann");
            }
            if(iter==1){
              int max_iteration, max_degree;
              Ifpack2::getParameter(test_params.sublist("Belos"), "Maximum Iterations", max_iteration);
              Ifpack2::getParameter(test_params.sublist("Ifpack2"), "neumann: degree", max_degree);
              max_iteration = static_cast<int>(2.0 * max_iteration /max_degree);
              test_params.sublist("Belos").set("Maximum Iterations", max_iteration);
            }
            
            //if(eig==2)
            //test_params.sublist("Ifpack2").set("eigen-analysis: type", "power method");
            //if(eig==2)
            //test_params.sublist("Ifpack2").set("neumann: degree", 2);
            //if(iter==2)
            //test_params.sublist("Ifpack2").set("neumann: degree", 2);

        } else if (prec == 2) {
            std::cout << mm_file <<" ==============================="<< std::endl;
            if(iter==2){
            //std::cout <<CYAN<< "----------- Neumann:ritz_updateSol_ritz -----------"<<RESET<< std::endl;
            std::cout <<CYAN<< "+++++++++++ poly_prec +++++++++++"<<RESET<< std::endl;
            }else{
            std::cout <<YELLOW<< "+++++++++++ poly_prec +++++++++++"<<RESET<< std::endl;
            }
            test_params.set("solver_type", "GmresPolySolMgr");
            test_params.set("Outer Solver", "Block Gmres");
            test_params.set("Ifpack2::Preconditioner", "ILUT");
            test_params.set("Preconditioner Side", "Right");
            test_params.sublist("Ifpack2").set("fact: ilut level-of-fill", 1.0);
            test_params.sublist("Ifpack2").set("fact: drop tolerance", 0.0);
            test_params.sublist("Ifpack2").set("fact: absolute threshold", 1.0e-8);
            if (test_params.sublist("Ifpack2").isParameter("neumann: degree")) {
              test_params.sublist("Ifpack2").remove("neumann: degree");
            }
            if(iter==1){
              int max_iteration, max_degree;
              Ifpack2::getParameter(test_params.sublist("Belos"), "Maximum Iterations", max_iteration);
              Ifpack2::getParameter(test_params.sublist("belosList"), "Maximum Degree", max_degree);
              //max_iteration = 2 * max_iteration /max_degree;
              max_iteration = static_cast<int>(2.0 * max_iteration /max_degree);
              test_params.sublist("Belos").set("Maximum Iterations", max_iteration);
            }
            
            //Teuchos::ParameterList outer_params=test_params.sublist("Belos");
            //test_params.set("Outer Solver Params",outer_params);
        }else if(prec == 4) {
            std::cout <<YELLOW<< "+++++++++++ Neumann:first_ritz +++++++++++"<<RESET<< std::endl;
            //    test_params.set("solver_type", "GmresPolySolMgr");
            //    test_params.set("Preconditioner Side", "None");
            test_params.set("solver_type", "Block Gmres");
            test_params.set("Ifpack2::Preconditioner", "NEUMANN");
            test_params.set("Preconditioner Side", "Right");

            if (!test_params.sublist("Ifpack2").isParameter("neumann: degree")|| iter==2) {
              //int degree;
              //Ifpack2::getParameter(test_params.sublist("Ifpack2"),"neumann: degree", degree); 
              test_params.sublist("Ifpack2").set("neumann: degree", 2);
            }
            //test_params.sublist("Ifpack2").set("neumann: degree", 1);
            test_params.sublist("Ifpack2").set("neumann: eigenvalue max iterations", 10);
            if (!test_params.sublist("Ifpack2").isParameter("eigen-analysis: type")) {
            test_params.sublist("Ifpack2").set("eigen-analysis: type", "power method");
            //test_params.sublist("Ifpack2").set("eigen-analysis: type", "ritz");
            }
            
            if (!test_params.sublist("Ifpack2").isParameter("neumann: algorithm")) {
            test_params.sublist("Ifpack2").set("neumann: algorithm", "textbook");
            }

        }else if (prec == 3) {       
            if(iter==2) continue;
          
            std::cout << mm_file <<" ==============================="<< std::endl;
            //std::cout << std::endl;
            //std::cout << "+++++++++++ ILU +++++++++++"<< std::endl;
            std::cout <<MAGENTA<< "============== ILUT =============="<<RESET<< std::endl;
            test_params.set("solver_type", "Block Gmres");
            test_params.set("Ifpack2::Preconditioner", "ILUT");
            test_params.set("Preconditioner Side", "Right");

            if (test_params.sublist("Ifpack2").isParameter("neumann: degree")) {
              test_params.sublist("Ifpack2").remove("neumann: degree");
            }
            test_params.sublist("Ifpack2").set("fact: ilut level-of-fill", 1.0);
            test_params.sublist("Ifpack2").set("fact: drop tolerance", 0.0);
            test_params.sublist("Ifpack2").set("fact: absolute threshold", 1.0e-8);
            
            //test_params.set("fact: absolute threshold", 1e-4);
        }else if (prec == 30) {
            std::cout << std::endl;
            std::cout <<CYAN<< "============== None =============="<<RESET<< std::endl;
            if (test_params.sublist("Ifpack2").isParameter("neumann: degree")) {
              test_params.sublist("Ifpack2").remove("neumann: degree");
            }
            test_params.set("solver_type", "Block Gmres");
            test_params.set("Ifpack2::Preconditioner", "RELAXATION");
            test_params.set("Preconditioner Side", "None");
         }else if (prec == 40) {
            std::cout << std::endl;
            std::cout <<CYAN<< "============== RELAXATION  =============="<<RESET<< std::endl;
            test_params.set("solver_type", "Block Gmres");
            test_params.set("Ifpack2::Preconditioner", "RELAXATION");
            test_params.set("Preconditioner Side", "Right");
            test_params.sublist("Ifpack2").set("relaxation: type", "Jacobi");

        }

    


    //} 
    //Ifpack2::getParameter(test_params, "mm_file", mm_file);
    
        RCP<BLinProb> problem =
          build_problem<Scalar,LO,GO,Node>(test_params, comm);
        //Teuchos::RCP<Tpetra::MultiVector<Scalar, L, GO, Node>> x0;
        Teuchos::RCP<TMV> x0,Av;
        x0 = problem->getLHS();
        RCP<TMV> R = rcp(new TMV(*problem->getRHS(), Teuchos::Copy));
        Teuchos::Array<Magnitude> normsR(R->getNumVectors());
        Teuchos::Array<Magnitude> normsB(R->getNumVectors());
        //x0->norm2(normsB);
        //R->norm2(normsR);
        //std::cout<<"x0_norm="<<normsB[0]<<std::endl;
        //std::cout<<"b_norm="<<normsR[0]<<std::endl;
        //->getVector(0)
 
        if(useStackedTimer)
        {
          stackedTimer = rcp(new StackedTimer("Belos-Ifpack2 Setup-Solve"));
          Teuchos::TimeMonitor::setStackedTimer(stackedTimer);
        }
 
        // defaullt convergence tol
        Magnitude tol = 1e-7;
        Teuchos::ParameterList& BelosParams = test_params.sublist("Belos");
        if (BelosParams.isParameter ("Convergence Tolerance")) {
          tol = BelosParams.get<Magnitude> ("Convergence Tolerance");
        } else {
          BelosParams.set<Magnitude> ("Convergence Tolerance", tol);
        }
 
        //Build the preconditioner (if one is enabled), and bind it to the problem
        std::string tifpack_precond("not specified");
        Ifpack2::getParameter (test_params, "Ifpack2::Preconditioner", tifpack_precond);
        std::string prec_side("Left");
        bool prec_info;
        Ifpack2::getParameter (test_params, "Preconditioner Side", prec_side);
        Ifpack2::getParameter (test_params, "Preconditioner Info", prec_info);
        if (tifpack_precond != "not specified") {
          RCP<TOP> precond;
          if (tifpack_precond == "RBILUK") {
            int blockSize = 0;
            Teuchos::ParameterList& prec_params = test_params.sublist("Ifpack2");
            Ifpack2::getParameter (prec_params, "fact: block size", blockSize);
            assert(blockSize >= 1);
            auto A_crs = Teuchos::rcp_dynamic_cast<const crs_matrix_type>(problem->getOperator());
            if (blockSize > 1) {
              auto crs_matrix_block_filled = Tpetra::fillLogicalBlocks(*A_crs, blockSize);
              auto A = Teuchos::rcp_const_cast<const block_crs_matrix_type>(Tpetra::convertToBlockCrsMatrix(*crs_matrix_block_filled, blockSize));
              precond = build_precond<Scalar,LO,GO,Node> (test_params, A, R, x0);
            }
            else {
              auto A = Teuchos::rcp_const_cast<const block_crs_matrix_type>(Tpetra::convertToBlockCrsMatrix(*A_crs, blockSize));
              precond = build_precond<Scalar,LO,GO,Node> (test_params, A, R, x0);
            }
          }
          else {
            auto A = Teuchos::rcp_dynamic_cast<const crs_matrix_type>(problem->getOperator());
            precond = build_precond<Scalar,LO,GO,Node> (test_params, A, R, x0);
          }
          if (prec_side == "Left")
            problem->setLeftPrec (precond);
          else if (prec_side == "Right")
            problem->setRightPrec (precond);
        }
 
        problem->setLHS(x0);
        problem->setProblem ();
 
        //The build_solver function is located in build_solver.hpp:
 
        RCP<BSolverMgr> solver = build_solver<Scalar,TMV,TOP>(test_params, problem,iter);
 
        Belos::ReturnType ret = solver->solve();
 
        *out << "   Converged in " << GREEN<< solver->getNumIters()<< RESET << " iterations." << std::endl;
 
        if (prec_info){
        *out << solver->description()<<std::endl;
        RCP<const TOP> prec = problem->getLeftPrec();
        if (prec !=Teuchos::null) {
          *out << "Preconditioner attributes:" << std::endl;
          prec->describe (*out, Teuchos::VERB_LOW);
        }
        prec = problem->getRightPrec();
        if (prec !=Teuchos::null) {
          *out << "Preconditioner attributes:" << std::endl;
          prec->describe (*out, Teuchos::VERB_LOW);
        }
        }
 
        //RCP<TMV> R = rcp(new TMV(*problem->getRHS(), Teuchos::Copy));
        problem->computeCurrResVec(&*R, &*problem->getLHS(), &*problem->getRHS());
        //Teuchos::Array<Magnitude> normsR(R->getNumVectors());
        //Teuchos::Array<Magnitude> normsB(R->getNumVectors());
        R->norm2(normsR);
        problem->getRHS()->norm2(normsB);
        
        //printf("normsB=,normsR=\n");
        //*out<<normsB<<", "<<normsR[0]<<std::endl;
 
        if (normsR.size() < 1) {
          throw std::runtime_error("ERROR: norms.size()==0 indicates R->getNumVectors()==0.");
        }
 
        //*out << "2-Norm of 0th RHS      vec: " << normsB[0] << std::endl;
        //*out << "2-Norm of 0th residual vec: " << normsR[0] << " -> " << normsR[0]/normsB[0] << std::endl;
        *out << "   Achieved tolerance: " << GREEN<<solver->achievedTol()<<", "<<normsR[0]/normsB[0]<<RESET << ", Requested tolerance: " << tol << std::endl;
        bool update_type(false);
        Ifpack2::getParameter(test_params, "updateSol", update_type);
        int off_set=0, degree = 2, iter_total;
        int iter_cnt=solver->getNumIters();
        double ratio;
        if(prec==1||prec==4){
          if(iter==1)
          Ifpack2::getParameter(test_params.sublist("Ifpack2"),"neumann: degree", degree); 
          if(update_type){off_set=11;}
          iter_total= degree*iter_cnt+off_set;
          ratio = exp(log(normsR[0]/normsB[0])/iter_total);
        }else if(prec==2){
          if(iter==1)
          Ifpack2::getParameter(test_params.sublist("belosList"), "Maximum Degree", degree);
          off_set=degree;
          iter_total= degree*iter_cnt+off_set;
          ratio = exp(log(normsR[0]/normsB[0])/iter_total);
        }else if(prec==3){
          iter_total= iter_cnt;
          ratio = exp(log(normsR[0]/normsB[0])/iter_total);
        }
        *out << "Matrix product count: " << GREEN<< iter_total << RESET << std::endl;
        *out << "Convegence Ratio:" << GREEN << ratio << RESET<<std::endl;
        normsR[0] /= normsB[0]; //=solver->achievedTol() when scaling=NormofRHS
 
        //If the xml file specified a number of iterations to expect, then we will
        //use that as a test pass/fail criteria.
 
        if (test_params.isParameter("expectNumIters")) {
          int expected_iters = 0;
          Ifpack2::getParameter(test_params, "expectNumIters", expected_iters);
          int actual_iters = solver->getNumIters();
          if (ret == Belos::Converged && actual_iters <= expected_iters && normsR[0] < tol) {
          }
          else {
            success = false;
            *out << "Actual iters(" << actual_iters << ") > expected number of iterations ("
                 << expected_iters <<"), or resid-norm(" << normsR[0] << ") >= " << tol << std::endl;
          }
        }
 
        timer.stop();
        auto timeT = timer.totalElapsedTime();
        double timeD = 100 * static_cast<double>(timeT);
        double perf=exp(log(ratio)*iter_total/timeD);
        *out <<RED<< "time: " << timeT << ", perf="<< perf <<RESET << std::endl;
        //*out <<RED<< "Performance: " << perf/(-log(normsR[0]))<<RESET << std::endl;

    }}//end, prec,eig
  }}//end iter,mm_file
    if(useStackedTimer)
    {
      stackedTimer->stopBaseTimer();
      StackedTimer::OutputOptions options;
      options.num_histogram=3;
      options.print_warnings = false;
      options.output_histogram = true;
      options.output_fraction=true;
      options.output_minmax = true;
      stackedTimer->report(std::cout, comm, options);
      auto xmlOut = stackedTimer->reportWatchrXML(problem_name, comm);
      if(comm->getRank() == 0)
      {
        if(xmlOut.length())
          std::cout << "\nAlso created Watchr performance report " << xmlOut << '\n';
      }
    }
  }
  
  TEUCHOS_STANDARD_CATCH_STATEMENTS(true, std::cerr, success)

  if (success) {
  //  *out << "End Result: TEST PASSED\n";
  }
  else {
  //  *out << "End Result: TEST FAILED\n";
  }

  return( success ? 0 : 1 );
}

